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Airborne  remote  sensing  for  archaeology  is  the  discipline  that  encompasses  the  study  of  archaeological 
remains  using  data  collected  from  an  airborne  platform  by  means  of  digital  or  film-based  aerial 
photography,  airborne  laser  scanning,  hyperspectral  imaging  etc.  So  far,  airborne  hyperspectral  scanning 
or  —  more  accurately  —  airborne  imaging  spectroscopy  (AIS)  has  occupied  only  a  very  small  niche  in  the 
field  of  archaeological  remote  sensing:  besides  reasons  of  cost,  the  common  archaeologically-insufficient 
ground-sampling  distance  can  be  considered  the  main  limiting  factor.  Moreover,  the  technical  processing 
of  these  data  sets  with  a  high  level  of  potential  redundancy  needs  specialized  software.  Typically, 
calculation  of  band  ratios  and  a  principal  component  analysis  are  applied.  As  a  result,  the  few  practical 
applications  of  archaeological  AIS  have  not  been  entirely  convincing  so  far.  The  aim  of  this  paper  is  to 
present  new  approaches  for  analysing  archaeological  AIS  data.  The  imagery  under  study  has  a  ground¬ 
sampling  distance  of  40  cm  and  covers  the  Roman  town  of  Carnuntum  (Austria).  Using  two  algorithms 
embedded  in  a  specifically  developed  MATLAB®  toolbox,  it  will  be  shown  how  the  extracted  archaeo¬ 
logical  information  can  be  enhanced  from  high-resolution  hyperspectral  images.  A  comparison  with 
simultaneously  acquired  vertical  photographs  will  indicate  the  specific  advantages  of  high-resolution  AIS 
data  and  the  gain  one  can  obtain  when  exploiting  its  potential  using  any  of  the  newly  presented 
methods. 

©  2014  The  Authors.  Published  by  Elsevier  Ltd.  This  is  an  open  access  article  under  the  CC  BY-NC-ND 

license  (http://creativecommons.Org/licenses/by-nc-nd/3.0/). 


1.  Introduction 

The  practice  of  aerial  photographic  reconnaissance  in  which  the 
archaeologist  acquires  predominantly  oblique  images  from  a 
manned,  low-flying  aeroplane  is  still  the  mostly  applied  archaeo¬ 
logical  remote  sensing  approach.  Although  this  survey  approach 
did  not  witness  major  changes  over  the  past  century  (Verhoeven, 
2009),  its  success  lies  mainly  in  its  straightforward  execution, 
while  its  capacity  to  cover  large  areas  turns  it  into  one  of  the  most 
cost-effective  methods  for  site  discovery.  The  non-invasive 
approach  also  yields  easily  interpretable  imagery  with  abundant 
spatial  detail  (Palmer,  2005;  Wilson,  2000).  The  latter  characteristic 


*  Corresponding  author.  Department  of  Prehistoric  and  Historical  Archaeology, 
University  of  Vienna,  Franz-Klein-Gasse  1,  1190  Wien,  Austria.  Tel.:  +43  1  4277 
40486:  fax:  +43  1  4277  9404. 

E-mail  address:  Michael. Doneus@univie.ac.at  (M.  Doneus). 


should  never  be  underestimated,  as  the  spatial  resolution  of  the 
acquired  datasets,  which  typically  is  between  5  cm  and  10  cm, 
enables  the  assessment  of  small,  but  often  abundant  features  such 
as  post  holes  and  pits.  It  can,  therefore,  be  considered  one  of  the 
main  reasons  that  this  type  of  aerial  survey  is  still  so  intensely 
practiced. 

Oblique  aerial  photography  does,  however,  exhibit  shortcom¬ 
ings  (discussed  below  in  more  detail),  mainly  due  to  its  limited 
spectral  resolving  power.  Therefore,  a  growing  number  of  archae¬ 
ologists  are  investigating  in  other  aerial  archaeology  techniques 
providing  images  with  increased  spectral  resolution.  Although 
airborne  multi-spectral  imaging  is  one  of  these  approaches  (e.g. 
Winterbottom  and  Dawson,  2005),  truly  increased  spectral  reso¬ 
lution  is  only  obtained  by  airborne  hyperspectral  scanning  (AHS) 
also  called  airborne  imaging  spectroscopy  (AIS)  (Aqdus  et  al„  2007, 
2008,  2012;  Barnes,  2003;  Bassani  et  al.,  2009a,  2009b;  Bennett 
et  al„  2012,  2013;  Cavalli  et  al.,  2007,  2012;  Cavalli,  2013;  Cavalli 
and  Pignatti,  2001;  Challis  et  al.,  2009;  Coren  et  al.,  2005; 
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Emmolo  et  al.,  2004;  Forte  et  al.,  2011;  Merola,  2005;  Merola  et  al., 
2007,  2008;  Pietrapertosa  et  al.,  2008;  Traviglia,  2005,  2006a, 

2006b,  2008;  White,  2003).  The  success  rate  in  terms  of  detection 
of  archaeological  subsurface  structures  is  varying,  however,  and 
less  successful  applications  seem  to  be  connected  with  the  lower 
spatial  resolution  of  the  acquired  datasets  (in  most  cases  ranging 
from  4  m  to  —  at  best  —  1  m).  Since  the  current  generation  of  A1S 
sensors  does  enable  the  acquisition  of  relatively  high  spatial  reso¬ 
lution  data  (i.e.  ground-sampling  distance  (GSD)  below  50  cm), 
they  seem  to  be  promising  for  archaeological  applications. 

The  high  spectral  resolution  of  AIS  sensors  is  generally  not  easily 
handled  in  standard  image  processing  software.  Indeed,  hyper- 
spectral  scans  result  in  huge  datasets  of  hundred  and  more  spectral 
bands,  which  are  moreover  to  a  high  degree  correlated  and  often 
affected  by  noise.  Processing  of  these  image  cubes,  i.e.  reducing  the 
data  to  only  a  few  layers  containing  a  high  degree  of  archaeologi- 
cally  relevant  information  is  therefore  an  important  part  of  the  AIS 
workflow. 

In  this  article,  we  will  present  two  new  approaches  that  have 
been  developed  and  tested  on  comparatively  high  spatial  resolu¬ 
tion  hyperspectral  datasets  (GSD  <50  cm).  They  are  implemented 
in  an  open  MATLAB®-based  toolbox  called  ARTC1S  (see  Atzberger 
et  al.,  2014)  and  make  use  of  the  spectral  features  stored  in  the 
AIS  data  cube.  Therefore,  a  basic  introduction  into  airborne 
hyperspectral  imaging  is  provided  in  the  next  section  (chapter  2), 
before  chapter  3  delves  deeper  into  the  two  new  processing  ap¬ 
proaches.  The  proposed  approaches  are  consequently  tested  on  a 
high-resolution  dataset  from  a  case-study  area  detailed  in  chapter 
4.  The  results  are  analysed  and  discussed  in  chapters  5  and  6. 

2.  Airborne  imaging  spectroscopy 

2.3.  Limitations  of  conventional  aerial  photography 

Despite  its  many  advantages  and  the  capability  to  capture  small 
landscape  details,  conventional  oblique  and  vertical  aerial  photog¬ 
raphy  is  also  characterised  by  some  major  spectral  shortcomings.  In 
oblique  aerial  photography,  the  survey  working  principles  are  gov¬ 
erned  by  the  human  visual  system  which  is  only  sensitive  to  the 
visible  electromagnetic  spectrum  between  400  and  700  nm.  Addi¬ 
tionally,  the  majority  of  the  aerial  footage  has  been  acquired  with 
photographic  media  that  were  sensitised  only  to  this  visible  radiation 
(e.g.  normal  colour  photography).  Although  very  striking  and 
revealing  images  have  been  obtained  in  this  way,  the  detection  of 
vegetation  marks  (and  to  a  certain  extent  soil  marks)  becomes 
impossible  in  less-optimal  circumstances,  as  the  slight  differentials 
of  height  and  colour  in  crops  might  exhibit  too  low  a  contrast  with 
the  surrounding  matrix  to  be  noticed  through  normal  (colour) 
photography  in  the  visible  spectrum.  Even  though  researchers  have 
been  experimenting  with  different  filters  and  film  emulsions 
(Crawshaw,  1995)  and  the  application  of  digital  cameras  has  enabled 
the  easy  application  of  beyond  visible  imaging  in  the  Near-InfraRed 
(NIR)  domain  (Verhoeven,  2012a),  visually  imperceptible  soil  and 
crop  disturbances  will  never  be  photographed.  As  an  answer  to  this 
‘observer-directed’  approach  (Palmer,  2005),  aerial  archaeologists 
have  also  used  photographs  generated  during  strictly  vertical  sorties 
(e.g.  Coleman,  2007;  Doneus,  1997;  Kennedy,  1996;  Mills,  2005; 
Palmer,  2007).  Although  even  NIR  coverage  becomes  possible  with 
such  a  vertical  approach,  both  the  standard  oblique  and  vertical 
photographic  strategies  still  capture  reflected  solar  radiation  in 
spectrally  broad  wavebands:  a  typical  colour  photograph  records  the 
visible  spectrum  in  three  approximately  100  nm  wide  channels  (Red, 
Green,  Blue),  whereas  NIR  radiation  is  usually  sampled  in  a  300  nm 
wide  band.  This  is,  however,  far  from  optimal,  because  particular 
diagnostic  spectral  features  are  often  only  a  few  nanometres  wide, 


which  makes  acquisition  of  data  with  a  high  spectral  resolution 
necessary  when  one  wants  to  assess  small  variations  in  —  for 
example  —  the  plants'  physiology.  Vertical  and  oblique  photography, 
even  when  executed  in  the  spectral  range  beyond  the  visible  light, 
thus  significantly  reduce  the  diagnostic  accuracy  of  vegetation 
investigation  (Hampton,  1974).  In  other  words,  they  hamper  the 
detection  of  archaeological  vegetation  marks  as  the  reflected  radia¬ 
tion  is  spectrally  undersampled  and  spectral  characteristics  that  are 
too  narrow  to  be  distinguished  get  masked  (Verhoeven,  2009). 

Spaceborne  data,  consistently  acquired  over  extended  areas  and 
often  in  invisible  wavebands,  have  been  used  in  a  variety  of 
archaeological  surveys  (Lasaponara  and  Masini,  2012).  They 
repeatedly  cover  large  parts  of  the  Earth's  surface,  are  often  easily 
available,  and  tackle  the  observer-directed  and  visible-radiation- 
limited  issues.  However,  the  data  are  less  (or  not  at  all)  suited  for 
the  discovery  and  detailed  recording  of  small  archaeological  fea¬ 
tures,  as  the  resolving  power  of  the  sensors  is  in  all  but  a  few  cases 
more  than  one  meter.  Moreover,  the  spectral  bands  of  older 
spaceborne  imagers  (i.e.  those  whose  products  are  freely  available) 
are  generally  too  broad  or  misplaced  spectrally  to  truly  detect  plant 
stress  (Carter,  1994).  Airborne  multi-  and  hyperspectral  sensors 
might  acquire  data  in  narrow  wavebands,  but  cost,  moderate 
temporal  resolution,  and  low  resolving  power  also  significantly 
hampered  their  frequent  use  in  archaeological  research  (Hanson, 
2008).  An  ideal  system  that  joins  the  best  of  these  approaches  by 
offering  the  cost-effectiveness  as  well  as  operating  and  post¬ 
processing  flexibility  of  the  oblique  reconnaissance  approach, 
while  also  allowing  a  total  coverage  in  narrow  visible  and  invisible 
spectral  wavebands,  does  not  exist  so  far.  However,  the  current 
generation  of  hyperspectral  sensors  does  enable  the  acquisition  of 
relatively  high  spatial  resolution  data  (i.e.  a  GSD  below  50  cm)  and 
due  to  the  increased  demand  in  non-archaeological  fields,  even  the 
acquisition  costs  are  steadily  decreasing.  However,  some  big  hur¬ 
dles  still  need  to  be  taken  before  these  high-resolution  AIS  data  can 
become  of  real  archaeological  interest. 

2.2.  AIS  —  what's  in  a  name? 

Airborne  hyperspectral  imaging  is  a  passive  remote  sensing 
technique  as  it  digitizes  the  earth's  upwelling  electromagnetic  ra¬ 
diation  (reflected  solar  radiation  or  thermal  radiation  emitted  by  the 
objects  themselves)  in  a  multitude  of  small  spectral  bands.  Just  as 
normal  photographs  consist  of  a  stack  of  three  spatially  co¬ 
registered  two-dimensional  images  of  which  each  represents  the 
reflected  radiation  in  a  broad  visible  band,  the  final  hyperspectral 
data  product  can  be  considered  an  extended  stack  of  narrow-band 
images.  Each  image  is  a  digitisation  of  the  reflected  and/or 
emitted  radiation  in  a  small  spectral  range  (typically  around  10  nm 
in  the  visible  and  NIR  range).  Commonly  an  airborne  hyperspectral 
flight  yields  tens  to  hundreds  of  these  narrowband  images,  which 
are  captured  in  spectrally  contiguous  bands.  The  final  product  can 
therefore  be  considered  a  three-dimensional  data  cube  (x,  y,  A)  in 
which  the  first  two  dimensions  are  the  spatial  dimensions,  whereas 
the  third  dimension  reflects  the  spectral  dimension  (Fig.  1).  Thus 
many  small  bands  of  electromagnetic  radiant  energy  are  captured 
per  pixel  location.  Just  as  a  pixel  of  a  common  digital  colour 
photograph  contains  three  samples  or  Digital  Numbers  (DNs)  at  the 
same  location  to  represent  the  amount  of  radiation  captured  in  the 
three  broad  spectral  bands,  a  hyperspectral  image  features  N  DNs,  in 
which  N  equals  the  amount  of  spectral  bands  that  are  sampled. 
Every  image  is  also  characterised  by  a  certain  bit  depth  (e.g.  integer 
10-bit  or  values  between  0  and  1023),  which  determines  the  reso¬ 
lution  by  which  the  at-sensor  radiance  L  can  be  mapped  onto  a 
discrete  set  of  digital  values.  From  these  DNs,  reflectance  (or  emis¬ 
sion)  can  be  calculated.  Through  a  combination  of  all  these  many 


86 


M.  Doneus  et  al.  /  Journal  of  Archaeological  Science  52  (2014)  84-96 


©Copyright Geert Verhoeven  (LBI  ArchPro) 


single  AIS 
raster  cell 


Fig.  1.  Visualization  of  an  AIS  data  cube. 


single  reflectance  values,  every  individual  pixel  of  the  final  image 
holds  the  complete  reflectance  or  emission  spectrum  (known  as 
spectral  signature)  of  the  material  that  was  sampled  at  that  specific 
location  (Fig.  1).  Since  this  spectral  signature  can  be  obtained  for 
every  pixel  in  the  image,  the  technique  is  also  (and  more  accurately) 
denoted  Airborne  Imaging  Spectroscopy  (AIS),  spectroscopy  being 
the  study  of  interaction  between  radiant  energy  and  matter. 

2.3.  AIS  information  extraction 

Because  imaging  spectroscopy  provides  a  spectrally  continuous 
sampling  of  the  objects  that  are  imaged,  it  allows  to  assess  specific 
and  subtle  vegetation  characteristics  (e.g.  stressed  versus  healthy 
plants)  and  soil  properties  (e.g.  mineral  composition)  to  a  much 
higher  extent  than  any  standard  photographic  or  multispectral 
method.  This  is  the  reason  why  archaeological  AIS  has  already 
several  years  ago  been  defined  as  a  domain  with  “huge  opportu¬ 
nities  for  the  efficient  detection,  mapping  and  management  of 
cultural  heritage”  (Traviglia,  2006b). 

Despite  the  better  spatial  and  spectral  resolutions,  working  with 
AIS  data  necessitates  expert  knowledge  and  specialised  software, 
which  is  again  a  limiting  factor  for  its  archaeological  use.  Assuming 
that  the  image  geometry  and  georeferencing  can  be  handled 
correctly,  one  needs  to  deal  with  the  data  quality  and  dimensional 
overload.  Indeed,  as  the  upwelling  electromagnetic  radiation  is 
recorded  in  small  spectral  bands,  the  signal  received  by  the  sensor 
is  quite  low  compared  to  sensor  noise  and  possible  atmospheric 
perturbations.  The  small  instantaneous  field  of  view  (IFOV)  of  the 


AIS  sensors  —  necessary  for  obtaining  a  sufficiently  high  spatial 
resolution  —  further  limits  the  useful  signal  incident  on  the  pho¬ 
todetectors.  To  deal  with  this  low  signal-to-noise  ratio  (SNR)  as  well 
as  the  huge  amount  of  available  information  that  is  not  directly 
accessible  by  the  human  visual  system,  data  mining  and  intelligent 
pre-  and  post-processing  approaches  are  necessary. 

To  date,  researchers  generally  apply  a  few  processing  techniques 
such  as  Principal  Component  Analysis  (PCA;  Aqdus  et  al.,  2008; 
Coren  et  al.,  2005;  Traviglia,  2006a)  and  various  vegetation 
indices  (Bennett  et  al.,  2012;  Cavalli  et  al.,  2007;  Emmolo  et  al., 
2004;  Traviglia,  2006b)  to  reduce  the  data  complexity  and  extract 
archaeologically  valid  information. 

In  essence,  PCA  is  a  mathematical  procedure  that  decreases  the 
data  redundancy  by  transforming  a  set  of  (possibly)  correlated 
variables  (original  image  data)  into  a  new  set  of  variables  that  are 
uncorrelated  —  the  latter  variables  denoted  the  Principal  Compo¬ 
nents  (PCs).  These  PCs  are  linear  combinations  of  the  original  var¬ 
iables  (Liu  and  Mason,  2009)  and  are  orthogonal  to  each  other.  The 
first  principal  component  (i.e.  the  variable  on  the  first  axis)  ac¬ 
counts  for  as  much  of  the  variability  in  the  original  data  as  possible, 
while  every  succeeding  principal  component  contains  decreasing 
amounts  of  variance  in  the  remaining  data.  The  first  PC  bands  thus 
holds  most  image  information,  while  the  higher  order  bands  are 
often  only  filled  with  noise,  although  they  still  might  display  very 
rare  features  (as  also  noted  by  Traviglia,  2006a).  In  remote  sensing, 
PCA  is  often  used  since  the  newly  calculated  component  images 
might  reveal  information  that  would  go  unnoticed  in  the  original 
data. 
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Vegetation  indices  (or  Vis)  are  arithmetic  combinations  of  DNs, 
radiance  or  reflectance  values  in  two  or  more  spectral  bands. 
Potentially,  this  combination  highlights  a  specific  feature  of  the 
vegetation  while  reducing  the  perturbing  effects  caused  by  solar 
geometry,  atmosphere,  topography,  viewing  angle  and  exposed 
soil.  There  exists  a  variety  of  Vis,  each  of  them  using  different 
datasets  and  optimised  for  specific  purposes  (e.g.  Bannari  et  al., 
1995).  As  a  result,  they  all  have  their  strengths  and  weaknesses  as 
well  as  operational  scales.  As  such,  not  all  these  Vis  are  of  archae¬ 
ological  benefit.  Although  disciplines  such  as  precision  farming 
might  require  specific  Vis  that  deliver  accurate  estimates  of  crop 
yield  and  leaf  chlorophyll  content,  archaeological  air-  and  space- 
borne  reconnaissance  rather  needs  reliable,  generally  applicable 
methods  that  can  be  extended  across  entire  landscapes,  at  different 
times  of  the  day,  and  during  various  growing  seasons  (rather  than 
very  specific  approaches  that  might  vary  from  field  to  field,  crop  to 
crop,  and  by  time  of  year  —  Verhoeven,  2012b). 

Although  both  PCA  and  Vis  have  proven  their  success  in  some 
cases,  more  advanced  processing  procedures  have  to  be  exploited 
now  that  AIS  data  sets  have  the  potential  to  even  disclose  very  small 
archaeological  features. 

3.  Case  study  and  data  acquisition 

To  test  the  potential  of  new  AIS  sensors  and  archaeologically- 
relevant  image  processing  approaches,  the  Ludwig  Boltzmann 
Institute  for  Archaeological  Prospection  and  Virtual  Archaeology 


(LBI  ArchPro  —  Trinks  et  al.,  2012)  had  set  up  a  strategy,  where  three 
Austrian  case  study  areas  —  Carnuntum,  Kreuttal  and  St.  Anna  — 
have  been  surveyed  repeatedly  between  2010  and  2012.  In  the 
following,  only  data  covering  the  area  of  Carnuntum,  located 
approximately  40  km  south-east  of  Vienna  (N  48°  6'  41",  E  16°  51' 
57"  —  WGS84),  was  used. 

Carnuntum  is  known  as  an  important  archaeological  area,  which 
is  situated  at  the  Roman  limes  on  the  southern  bank  of  the  Danube 
River  (Fig.  2).  It  was  an  important  focal  point  throughout  the  first 
four  centuries  AD  and  as  such  the  capital  of  the  former  Roman 
province  Pannonia  superior  (Jobst  et  al.,  1983).  The  area  is  charac¬ 
terized  by  the  river  Danube  and  its  gravel  terraces.  Buried  under  an 
intensively  cultivated  calcic  chernozem  (Neubauer  et  al.,  2002),  the 
archaeological  remains  covering  an  area  of  roughly 
6.5  square  kilometres  usually  produce  well  visible  vegetation 
marks  during  the  month  of  June.  During  the  last  decades,  Car¬ 
nuntum  has  been  a  focus  for  non-destructive  archaeological  pro¬ 
spection  (Doneus  et  al.,  2013;  Neubauer  et  al.,  2012,  2014;  Saey 
et  al.,  2013;  Seren  et  al.,  2007).  Specifically,  since  2010  it  has 
become  a  case-study  area  of  the  LBI  ArchPro.  As  a  result,  a  complete 
interpretation  of  the  aerial  archaeological  evidence  from  50  years  of 
reconnaissance  is  available  (Doneus  et  al.,  2013). 

Altogether,  four  scanning  campaigns  were  realized  by  the  Aus¬ 
trian  LBI  partner  Airborne  Technologies.  All  of  the  data  were  ac¬ 
quired  in  high  spatial  and  spectral  resolution  with  GSDs  between 
40  cm  and  60  cm  in  65—105  spectral  bands  between  400  nm  and 
1000  nm  and  digitised  as  12-bit  integer  DNs  (Table  1). 


Fig.  2.  The  case  study  area  of  Carnuntum.  1:  Legionary  fortress;  red  rectangle:  sample  area  described  in  the  discussion  below  (see  also  Figs.  5  and  7).  (For  interpretation  of  the 
references  to  colour  in  this  figure  legend,  the  reader  is  referred  to  the  web  version  of  this  article.) 
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Table  1 

Flight  parameters  of  the  AIS  data  acquisition  campaigns  over  the  case  study  area  of  Carnuntum,  Lower  Austria.  (FWHM  means  Full  Width  at  Half  Maximum). 


Case  study  area 

Carnuntum 

Carnuntum 

Carnuntum 

Carnuntum 

Date  of  data  acquisition 

5  June  2010 

26  May  2011 

11  November  2011 

18  June  2012 

Purpose  of  scan 

Archaeology 

Archaeology 

Archaeology 

Archaeology 

Imaging  spectrometer 

AisaEAGLET  (SPECIM) 

AisaEAGLE  (SPECIM) 

AisaEAGLE  (SPECIM) 

AisaEAGLE  (SPECIM) 

Number  of  bands 

105 

65 

65 

65 

Spectral  sampling  interval 

53—5.9  nm 

8.6— 9.5  nm 

8.6— 9.5  nm 

8.6— 9.5  nm 

Spectral  resolution 

3.3  nm  at  FWHM 

3.3  nm  at  FWHM 

3.3  nm  at  FWHM 

3.3  nm  at  FWHM 

Spectral  range 

400-1001  nm 

403—995  nm 

403-995  nm 

403—995  nm 

Binning 

8 

8 

8 

8 

Flying  height  above  ground 

550/700  m 

600  m 

600  m 

600  m 

Ground  Sampling  Distance 

0.4  m/0.6  m 

0.4  m 

0.4  m 

0.4  m 

Speed  of  aircraft 

50  m/sec 

50  m/sec 

50  m/sec 

50  m/sec 

Digitisation 

12-bit 

12-bit 

12-bit 

12-bit 

Simultaneously  with  the  hyperspectral  data  acquisition  (using 
an  AISA  Eaglet  and  AISA  Eagle  device),  a  full  waveform  laser  scan  (at 
least  10  points  per  m2)  using  a  RIEGL  LMS-Q680i  and  vertical 
photographs  (IGl-DigiCam  H39  with  a  50  mm  lens)  were  acquired. 
The  laser  scan  resulted  in  a  high  quality  digital  surface  model 
(DSM)t  which  was  used  for  the  orthorectification  of  both  the 
hyperspectral  scan  and  the  vertical  aerial  photographs. 

Both  atmospheric  correction  and  accurate  georeferencing  were 
necessary  to  combine  the  resulting  multi  temporal  hyperspectral 
data  and  to  integrate  them  with  other  prospection  techniques 
applied  by  the  LBI  Archpro  (mainly  aerial  photography,  magnetics 
and  ground-penetrating  radar).  This  was  achieved  using  the  ap¬ 
plications  ATCOR  for  radiometric  correction  (only  for  AisaEAGLE 
data)  and  PARGE®  (http://www.rese.ch)  for  accurate  georeferenc¬ 
ing  with  sub-pixel  accuracy. 

4.  New  ways  of  AIS  information  extraction  using  the  ARCTIS 
toolbox 

To  analyse  the  corrected  and  orthorectified  AIS  data  blocks,  a 
MATLAB®  toolbox  called  ARCTIS  (ARChaeological  Toolbox  for  Im¬ 
aging  Spectroscopy)  with  accompanying  user-friendly  Graphical 
User  Interface  (GUI)  was  developed  (Atzberger  et  al.,  2014; 
Atzberger  and  Wess,  2012).  The  toolbox  allows  the  image  analyst 
(not  necessarily  a  specialist  in  remote  sensing  or  imaging  spec¬ 
troscopy)  to  get  a  maximum  of  archaeological  information  out  of 
the  recorded  3D  data  cube.  As  the  main  focus  for  the  toolbox  is 
archaeological  prospection,  the  aim  was  to  visualize  the  data  (or 
subsets  of  it)  and  extract  possibly  occurring  vegetation  and  soil 
marks. 

The  toolbox  allows  to  analyse  small  datasets  providing  all  of  the 
functions  which  have  been  identified  in  the  past  as  useful  for 
archaeological  information  extraction.  In  addition,  several  new 
functionalities  were  developed  and  integrated.  In  a  first  step  the 
geometrically  and  radiometrically  corrected  data  are  imported  and 
some  further  pre-processing  is  performed.  To  deal  with  the  noise  in 
the  individual  spectral  bands  (which  is  even  more  present  when  no 
radiometric  correction  is  applied  beforehand),  a  very  powerful  and 
commercially  unavailable  spectral  smoothing  filter  (based  on  the 
Whittaker  smoother  —  Atzberger  and  Eilers,  2011a,  2011b)  was 
implemented  in  ARCTIS.  After  its  executing,  a  very  great  deal  of 
image  noise  is  removed  from  all  individual  spectral  bands  without 
risking  the  removal  of  important  spectral  features  or  inducing  new 
spectral  artefacts. 

For  the  subsequent  information  extraction  (from  now  on  simply 
denoted  as  processing),  several  tools  have  been  compiled;  PCA,  a 
collection  of  the  most  important  vegetation  indices,  a  similarity 
function,  edge-extraction  techniques,  and  classification  by  k-means 
clustering.  Due  to  the  highly  correlated  nature  of  the  spectral 


information  in  the  individual  bands  of  an  AIS  dataset  and  the  omni¬ 
present  image  noise,  image  processing  techniques  such  as  PCA  or 
Minimum  Noise  Fraction  (MNF)  are  often  used  to  calculate  a 
spectral  subset  of  the  initial  data  cube  by  rejecting  the  noisy 
components.  The  two  new  methods  presented  in  the  following 
paragraphs  —  red  edge  inflection  point  (REIP)  and  distribution 
fitting  —  take,  however,  the  complete  AIS  dataset  into  account:  REIP 
describes  the  spectral  shape  of  the  reflectance  curve  whereas  the 
distribution  fitting  parameterizes  the  shape  of  the  frequency  dis¬ 
tribution  of  the  reflectance  values. 

4.1  Red  edge  inflection  point 

In  the  first  method,  called  the  Red  Edge  Inflection  Point  (REIP), 
specific  attention  is  given  to  the  Red  Edge.  This  Red  Edge  is  a 
transition  zone  around  700  nm  in  the  visible  far-Red  to  NIR  tran¬ 
sition  spectrum,  for  healthy  vegetation  characterized  by  a  steep 
slope  of  the  spectral  reflectance  curve  (Fig.  3A).  This  very  abrupt 
reflectance  change  results  from  the  fact  that  the  absorption  by 
chlorophyll  pigments  in  the  far-Red  spectrum  (around 
670—680  nm)  is  high  and  the  NIR  reflection  by  the  spongy  meso- 
phyll  is  high,  giving  rise  to  one  of  the  most  extreme  slopes  to  be 
found  in  reflectance  spectra  of  natural  materials.  When  vegetation 
gets  stressed,  this  Red  Edge  shifts  toward  shorter  wavelengths  due 
to  the  stress  related  increase  of  reflectance  at  wavelengths  near 
700  nm  (Fig.  3A).  This  “Blue  shift”  of  the  Red  Edge  is  currently 
accepted  as  the  most  consistent  stress  response  of  plants  to 
different  stressors  (Boochs  et  al.,  1990;  Gitelson  et  al.,  1996;  Horler 
et  al.,  1983).  Previous  research  has  already  proven  the  archaeo¬ 
logical  usefulness  of  the  red  edge  spectral  band  when  used  in  a 
simple  vegetation  index  (Bennett  et  al„  2012;  Verhoeven,  2009; 
Verhoeven  and  Doneus,  2011;  Cho  et  al.,  2008;  Darvishzadeh 
et  al.,  2009).  Since  AIS  sensors  often  sample  the  Red  Edge  in  a 
few  bands,  it  is  one  of  the  first  spectral  regions  one  should  focus  on 
when  looking  for  vegetation  stress  in  the  complete  optical  region 
(Carter  et  al.,  2002). 

The  REIP  algorithm  in  ARCTIS  will  use  the  complete  red  edge 
spectrum  to  look  for  its  inflection  point  (as  is  marked  by  a  red  cross 
in  the  spectral  curve  in  Fig.  3B).  To  this  end,  the  user  first  defines  the 
region  he  wants  the  algorithm  to  focus  on.  For  the  REIP,  often  the 
680  nm— 760  nm  region  is  chosen.  However,  this  freedom  of  choice 
indicates  that  one  is  not  limited  to  the  red  edge.  Although  not 
denoted  REIP  anymore,  inflection  points  might  also  be  estimated  in 
the  short-wavelength  infrared  or  the  green  waveband.  This  inflec¬ 
tion  point  is  the  location  of  the  highest  gradient  in  the  spectral 
profile  of  the  predetermined  spectral  range. 

In  ARCTIS,  some  pre-processing  steps  can  be  performed  to 
obtain  more  accurate  results  (Fig.  3B)  (Atzberger  et  al.,  2014).  To 
start,  the  Whittaker  smoother  can  be  applied  during  the  call  of  the 
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Fig.  3.  (A)  Typical  reflectance  responses  of  a  Radiata  pine  leaf.  The  discontinuous  green  curve  represents  the  mean  reflectance  for  healthy,  chlorophyll-rich  leaves,  while  the  upper 
orange  curve  indicates  the  mean  reflectance  for  a  stress  situation  (nitrogen  deficiency).  The  grey  curve  represents  the  standard  deviation  of  reflectance  calculated  from  reflectance 
spectra  of  maize  leaves  with  varying  chlorophyll  content  (Verhoeven  and  Doneus,  2011,  Fig.  1  A).  (B)  The  Red  Edge  Inflection  Point  GUI  as  it  is  currently  implemented  in  the  ARCTIS 
toolbox  (see  text  for  more  details  on  the  available  functionality).  The  red  cross  marks  the  position  of  the  inflection  point.  (For  interpretation  of  the  references  to  colour  in  this  figure 
legend,  the  reader  is  referred  to  the  web  version  of  this  article.) 


REIP  algorithm,  where  each  A1S  pixel  is  smoothed  independently 
along  the  z-axis  (i.e.  spectral  dimension)  according  to  a  user- 
defined  smoothing  parameter  lambda  (7).  Additionally,  the 
smoothed  input  AIS  data  can  be  spectrally  interpolated  to  a  user- 
defined  number  of  fictional  bands  between  the  original  spectral 
bands  leading  to  a  finer  spectral  resolution  of  the  detected  REIP 
location  (Fig.  3B).  Whatever  pre-processing  steps  were  chosen,  the 
REIP  is  finally  calculated  on  a  pixel-by-pixel  basis,  where  for  each 
pixel  the  inflection  point  is  located,  i.e.  the  steepest  slope  of  the 
smoothed  and  interpolated  spectral  curve  within  the  predefined 
spectral  region.  From  each  inflection  point,  three  parameters  can  be 
derived,  which  can  be  compiled  into  three  new  image  layers: 

•  Layer  1:  band  or  wavelength  of  the  REIP  location; 

•  Layer  2:  value  of  the  slope  of  the  reflectance  curve  at  the  REIP 
position; 

•  Layer  3:  reflectance  value  or  DN  at  REIP. 


4.2.  Distribution  fitting 

A  second  technique  embedded  in  ARCTIS  is  called  distribution 
fitting,  for  which  an  algorithm  tries  to  fit  a  user-defined  function  to 
the  frequency  distribution  of  a  pixel's  spectral  values  (Atzberger 
et  al.,  2014).  Fig.  4  depicts  this  process.  First,  the  algorithm  re¬ 
trieves  all  the  spectral  data  from  one  hyperspectral  pixel.  The  latter 
consists  of  various  reflectance  values:  one  for  each  of  the  acquired 
spectral  bands.  Since  these  data  are  univariate  (i.e.  they  provide 
information  on  one  variable,  being  reflectance)  a  frequency  histo¬ 
gram  can  be  calculated.  Fig.  4  shows  the  histograms  and  spectral 
signatures  of  two  completely  different  AIS  pixels:  one  pixel  stores 
the  reflectance  values  from  a  healthy  green,  chlorophyll-rich  leaf, 
while  the  other  pixel  digitised  the  reflectance  of  a  stressed,  yellow 
leaf.  From  the  spectral  reflectance  curve,  it  is  obvious  that  the 
stressed  leaf  suffered  from  a  severe  loss  of  Chlorophyll  pigments. 
First  of  all,  the  Green  reflectance  peak  around  550  nm  is  broadened 
towards  longer  wavelengths  and  causes  the  leaf  to  appear 
yellowish  (a  phenomenon  called  chlorosis  —  Adams  et  al.,  1999).  At 
the  same  time,  the  Red  Edge  shifts  toward  shorter  wavelengths  due 
to  the  stress  related  increase  of  reflectance  at  wavelengths  near 
700  nm  (i.e.  the  Blue  shift).  The  shape  of  the  spectral  curves  is  also 


directly  influencing  their  corresponding  histogram.  Consider  for 
example  the  reflectance  in  the  25%— 30%  interval.  From  the  spectral 
curve,  one  can  infer  that  the  healthy  green  leaf  has  only  one  band 
(band  number  53  in  the  red  edge  region)  with  a  reflectance  value  in 
this  range,  which  leads  to  a  single  count  in  this  bin.  This  contrasts 
highly  with  the  chlorotic  leaf  pixel,  which  consists  of  almost  twenty 
bands  (i.e.  bands  30—50)  with  a  reflectance  value  in  this  range.  As  a 
result,  the  histogram  bin  for  this  reflectance  range  has  a  height  of 
eighteen  counts.  This  relationship  between  the  reflectance  curve 
and  its  histogram  now  clarifies  the  more  homogenous  frequency 
distribution  in  the  chlorotic  histogram  with  only  one  large  NIR 
peak,  whereas  the  healthy  leaf  histogram  is  dominated  by  two  large 
peaks  that  are  related  to  the  low  reflectance  in  the  visible  range  and 
the  high  NIR  reflectance  respectively. 

After  building  the  frequency  distribution  histogram  for  one  AIS 
pixel,  a  predefined  Probability  Density  Function  (PDF)/is  fitted  to  it. 
In  probability  theory,  various  PDFs  are  used  to  mathematically 
describe  the  specific  distribution  of  a  random  variable  (the  reflec¬ 
tance  in  this  case).  Since  random  variables  can  be  continuous  or 
discrete,  both  types  of  PDFs  exist  as  well.  In  the  continuous  cate¬ 
gory,  many  different  PDFs  reside:  normal  or  Gaussian,  Rayleigh, 
exponential,  Weibull,  beta,  gamma,  chi,  Student’s  t-distribution  etc. 
The  Poisson  and  Boltzmann  density  functions  are  both  examples  of 
the  discrete  class.  Any  of  these  PDFs  is  defined  by  its  parameters. 
The  normal  density  function  describing  the  normal  distribution,  for 
instance,  is  parameterized  in  terms  of  a  mean  parameter  n  and  a 
standard  deviation  a  (Moore  and  McCabe,  2003).  In  ARCTIS,  a 
number  of  discrete  and  continuous  theoretical  PDFs  reside  and  all 
can  be  parameterized. 

Thus,  after  generating  the  histogram  from  the  acquired  signal 
(i.e.  the  pixel's  values),  the  algorithm  will  parameterise  the  chosen 
PDF  so  that  it  best  characterizes  this  pixel-specific  histogram.  The 
values  of  the  parameters  describing  this  PDF  will  then  make  up  the 
DNs  of  the  newly  generated  image;  the  number  of  image  layers  will 
equal  the  number  of  parameters  of  the  PDF.  Fig.  4  shows  the  results 
one  gets  with  two  different  PDFs:  the  normal  and  the  gamma 
continuous  density  functions.  Since  both  PDFs  are  mathematically 
different,  they  will  obviously  yield  dissimilar  results  for  the  same 
histogram.  Besides,  different  histograms  will  yield  dissimilar  pa¬ 
rameters  for  the  same  PDF.  So,  in  the  case  of  the  green  leaf,  a  new 
two-band  image  will  be  created  with  pixel  values  31  and  21  when 
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Fig.  4.  The  principle  of  a  pixel-wise  distribution  fitting  (more  information  provided  in  the  text). 


choosing  the  normal  PDF.  This  process  is  then  repeated  for  each 
pixel  of  the  A1S  data  block.  When  the  neighbouring  pixel  of  the 
original  AIS  image  holds  for  instance  the  reflectance  spectrum  of 
the  chlorotic  leaf,  the  new  pixel  values  will  become  33  and  14. 

The  result  of  this  very  computing  intensive  process  is  a  two- 
band  image,  in  which  the  difference  between  stressed  and 
healthy  vegetation  is  visualised  in  a  completely  new  manner.  The 
exact  amount  of  image  channels  that  are  computed  completely 
depends  on  the  PDF  selected.  Obviously,  in  the  case  one  would  pick 
a  PDF  that  is  parameterised  by  three  values,  a  three-band  image  is 
produced.  Initial  promising  results  obtained  with  this  new  tech¬ 
nique  were  briefly  reported  by  Verhoeven  et  al.  (2013),  but  here  the 
coverage  will  be  more  in  depth. 

After  testing  a  handful  of  functions,  the  gamma  and  normal 
continuous  density  functions  yielded  the  best  results  on  our 
dataset.  This  is  due  to  the  fact  that  the  gamma  density  function 
can  cope  very  well  with  the  right  skewed  histograms  from 
healthy  vegetation,  which  makes  its  parameterisation  completely 
different  when  dealing  with  stressed  vegetation  histograms  (see 
Fig.  4).  Although  the  gamma  density  function  is  often  parame¬ 
terized  with  a  shape  parameter  k  and  scale  parameter  8,  ARCTIS 
uses  the  parameterization  embedded  in  MATLAB®,  which  is 
based  on  a  shape  parameter  a  (=  k)  and  a  rate  parameter  b  (=  1/8) 
(The  MathWorks,  2013).  The  Normal  density  function  is  param¬ 
eterized  in  terms  of  a  mean  parameter  p  and  a  standard  deviation 


a  (Moore  and  McCabe,  2003).  Since  MATLAB  produces  two 
additional  values,  being  the  lower  and  upper  bound  of  the  con¬ 
fidence  interval  for  p  (The  MathWorks,  2014),  the  results  is  a  four 
band  image.  One  does,  however,  not  have  to  use  all  four  gener¬ 
ated  bands  and  can  only  focus  on  one  band  (i.e.  a  greyscale  im¬ 
age)  or  on  two  or  three  bands  (i.e.  a  false-colour  image).  Since  the 
mean  p  of  the  normal  PDF  is  often  not  very  dissimilar  for  stressed 
versus  healthy  vegetation  fittings  (see  Fig.  4),  its  standard  devi¬ 
ation  <7  often  displays  the  contrast  best  (see  also  next  section  and 
Fig.  5C). 

5.  Results 

For  the  following  discussion,  a  sample  area  was  chosen  within 
the  former  canabae  legionis  south-west  of  the  Roman  legionary 
fortress  (see  red  rectangle  in  Fig.  2).  It  contains  several  fields  with 
different  crops. 

Fig.  5  clearly  indicates  the  remarkable  enhancement  one  can 
obtain  through  REIP  processing  and  distribution  fitting.  Data 
acquisition  was  on  May  26  2011,  which  was  due  to  wet  conditions 
during  April  at  the  beginning  of  the  “vegetation-mark  season”, 
when  differences  in  height  and  colour  of  the  vegetation  started  to 
emerge.  According  to  Metop/ASCAT  soil  moisture  data  viewer 
(ASCAT,  2014)  the  degree  of  saturation  of  the  soil  moisture  content 
of  the  upper  soil  layer  (<2  cm)  was  12.55%. 
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Fig.  5.  Data  from  Carnuntum,  Lower  Austria,  acquired  on  May  26  2011.  (A)  Conventional  orthorectified  aerial  image,  acquired  in  the  visible  spectrum.  GSD  of  0.1  m;  enhanced  using 
contrast  limited  adaptive  histogram  equalization  (CLAHE);  (B)  false  colour  composite  created  by  means  of  the  RE1P  algorithm  {R  =  band  1  (wavelength),  G  =  band  2  (slope),  B  =  band 
3  (reflectance  value)),  (C)  rate  parameter  b  of  the  gamma  distribution  fitting;  (D)  normal  distribution  fitting  (  R  =  NONE,  G  =  band  2  (rr),  B  =  band  1  (p)).  GSDs  of  (B),  (C)  and  (D): 
0.4  m.  Figures  B,  C,  and  D  were  subject  to  the  same  histogram  stretch  by  means  of  standard  deviation.  (For  interpretation  of  the  references  to  colour  in  this  figure  legend,  the  reader 
is  referred  to  the  web  version  of  this  article.) 


To  improve  the  visibility  of  archaeological  traces  captured  in  the 
conventional  aerial  photograph,  a  contrast  limited  adaptive  histo¬ 
gram  equalization  ( CLAHE  —  Zuiderveld,  1994)  was  applied 
(Fig.  5A).  Still  the  photograph  displays  only  very  few  marks,  which 
are  almost  entirely  located  in  the  narrow  dark-green  field  extend¬ 
ing  from  the  upper  left  corner  towards  the  middle  of  the  lower 
edge.  They  can  be  described  as  negative  vegetation  marks,  which 
are  a  consequence  of  (in  this  case)  buried  walls  of  Roman  buildings. 
In  the  large  field  to  the  right,  a  Roman  road  is  revealed  by  a  linear 
arrangement  of  yellow  (i.e.  already  stressed)  crops  (Fig.  6A). 

The  false  colour  composite  created  by  means  of  the  RE1P  algo¬ 
rithm  (R  =  wavelength,  G  =  slope,  B  =  reflectance  value)  is  dis¬ 
played  in  Fig.  5B.  It  enhances  the  differences  in  plant  physiology 
and  consequently  allows  the  archaeological  features  to  be  seen 
more  clearly.  Especially  the  above  mentioned  negative  vegetation 
marks  in  the  narrow  field  can  now  be  read  in  detail,  while  further 
archaeologically  induced  traces  become  visible  in  the  larger  field  to 
the  right.  Their  contrast  is  fainter,  but  they  can  be  identified  as 
remains  of  buried  buildings  located  on  both  sides  of  the  Roman 
road,  which  is  now  clearly  visible  in  dark  red  (Fig.  6B). 

The  results  of  fitted  gamma  (Fig.  5C)  and  normal  (Fig.  5D)  PDFs 
perform  almost  equally  well.  For  the  gamma  distribution  fitting,  the 
histogram-stretched  version  of  the  rate  parameter  b  gave  the  best 
results,  while  a  combination  of  the  standard  deviation  a  (band  2) 
and  the  upper  bound  of  the  confidence  interval  for  ji  (band  4) 
(R  =  band  4,  G  =  band  2,  B  =  band  2)  turned  out  to  be  the  best 
combination  to  identify  vegetation  marks  in  this  scene.  The 
explanation  for  this  can  be  deduced  from  the  spectral  signatures  in 
Fig.  4,  as  they  are  very  similar  to  the  spectral  reflectance  curves  of 
the  chlorotic  and  green  vegetation  in  these  fields.  Both  fittings 


clearly  indicate  negative  vegetation  marks  in  darker  tones. 
Although  vegetation  marks  seem  to  be  slightly  enhanced  in  the  PDF 
images,  in  terms  of  archaeological  features  there  is  virtually  no 
difference  to  the  REIP  image  (Fig.  6C  and  D  compared  to  Fig.  6B). 
Slight  variations  do  occur,  where  for  each  visualization  other  details 
have  better  visibility. 

The  situation  is,  however,  different,  when  the  season  gets 
warmer  and  drier.  Fig.  7  shows  the  slightly  larger  and  shifted  area  of 
Fig.  5  as  recorded  on  June  18,  2012,  after  a  rather  dry  spring  (ZAMG, 
2014).  The  relative  dry  conditions  are  also  reflected  in  the  soil 
moisture  data  (ASCAT,  2014),  which  were  relatively  low  during 
April  and  May  and  display  around  14%  saturation  of  the  soil  mois¬ 
ture  content  of  the  upper  soil  layer  (<2  cm)  at  the  date  of  data 
acquisition. 

Both  visualizations  of  the  REIP  image  and  gamma  PDF  have  the 
same  band  combinations  as  in  Fig.  5.  For  the  normal  PDF  the 
standard  deviation  o  was  chosen  as  the  best  visualization  in  this 
case  (Fig.  7D).  In  the  lower  left  part  of  the  image,  the  Roman  road 
from  Fig.  5  is  again  visible,  but  due  to  a  different  crop  (maize),  no 
traces  of  the  buildings  can  be  seen  in  any  of  the  visualizations.  The 
two  fields  (annotations  1  and  2  in  Fig.  7A)  to  the  west  do  show 
detailed  negative  vegetation  marks.  Due  to  a  soil  moisture  deficit, 
the  summer  crop  in  field  1  shows  a  yellow  on  green  contrast  be¬ 
tween  stressed  and  healthy  plants,  which  is  extremely  well  visible 
in  the  conventional  photograph  (Fig.  7A).  This  situation  is  more  or 
less  similar  to  the  one  displayed  in  Fig.  5,  although  slightly  more 
pronounced  in  this  case.  Field  2  is  grown  with  winter  crop,  which 
has  already  entirely  ripened.  Negative  vegetation  marks  show  in  a 
darker  brown  colour,  which  makes  a  (less  well  visible)  contrast  to 
the  surrounding  yellow  plants. 
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Fig.  6.  Mapping  of  visualizations  from  Fig.  5:  (A)  orthophotograph,  (B)  RE1P,  (C)  gamma  distribution  fitting,  (D)  normal  distribution  fitting. 


In  this  situation,  when  vegetation  marks  become  clearly  evident 
in  the  visible  part  of  the  electromagnetic  spectrum,  the  false  colour 
image  produced  by  the  REIP  algorithm  cannot  add  or  enhance 
archaeological  information:  while  vegetation  marks  in  field  1  are 
slightly  less  well  visible  (which  is  mainly  due  to  the  four  times 
lower  spatial  resolution),  calculation  of  the  REIP  is  not  successful  at 
all  in  the  already  ripened  vegetation  of  field  2  (Fig.  SB).  This  con¬ 
firms  what  was  already  indicated  in  previous  research:  taking  the 
Red  edge  into  account  is  only  very  beneficial  when  one  wants  to 
identify  and  enhance  weak  vegetation  stresses  that  are  lost  when 
taking  only  the  broad  visible  spectrum  into  account  (Verhoeven, 
2009;  Verhoeven  and  Doneus,  2011).  When  vegetation  is  imaged 
at  the  very  end  of  the  season  and  all  plants  are  senescing,  the  red 
edge  has  completely  disappeared.  As  a  result,  the  REIP  approach 
uses  its  edge  over  the  visible  and  pure  NIR  frames.  Both  gamma 
(Fig.  8C)  and  normal  (Fig.  8D)  distribution  function  show  better 
results  than  the  REIP  algorithm,  but  also  cannot  outperform  the 
conventional  photograph  during  this  part  of  the  growing  season. 

6.  Discussion 

These  two  examples  show  that  airborne  hyperspectral  scanning 
can  reveal  detailed  archaeological  information,  if  a  high-resolution 
data  acquisition  with  a  GSD  below  1  m  is  available.  At  times  when 
the  vegetation  mark  season  is  at  its  beginning,  or  in  a  rather  wet 
growing  season  (as  was  the  case  in  2011  —  see  Fig.  5)  hyperspectral 
data  may  even  outperform  conventional  aerial  photography.  A  high 
spatial  resolution  will,  however,  increase  the  noise  in  the  individual 
spectral  bands,  which  needs  sophisticated  preparation  of  the  data 
cube.  In  our  case  the  Whittaker  smoother,  a  very  powerful  spectral 


smoothing  filter,  was  applied  (for  examples  with  and  without 
smoothing,  please  consult  Atzberger  et  al.,  2014). 

Both  introduced  methods  of  information  extraction  —  REIP  and 
distribution  fitting  —  are  useful  add-ons  to  PCA  and  vegetation 
indices,  which  are  the  currently  used  techniques  for  reduction  of 
data  complexity  and  extraction  of  archaeological  information.  For 
the  sake  of  this  paper,  mappings  of  all  visualizations  (Figs.  6  and  8) 
were  added  to  demonstrate  the  archaeological  information  content 
of  each  method.  A  detailed  assessment  of  the  qualitative  and 
quantitative  performance  of  all  these  different  methods  as  well  as  a 
comparison  with  other  techniques  (e.g.  vegetation  indices  or  PCA) 
will  be  however  subject  of  another  paper.  Also  further  systematic 
testing  in  different  environments  and  different  archaeological  site 
types  will  be  necessary. 

The  results  of  the  REIP  and  distribution  fittings  (as  is  also  the 
case  with  PCA)  are  new  bands,  which  are  displayed  by  the  imaging 
software  (in  our  case  ArcGIS)  as  false-colour  colour  or  greyscale 
images.  Depending  on  the  situation  of  the  data  acquisition  (land 
cover,  date  and  time  of  data  acquisition  etc.)  and  the  purpose  and 
method  of  interpretation,  different  assignments  of  the  resulting 
bands  with  one  of  the  three  image-channels  will  be  necessary.  This 
can  be  seen  in  this  paper  with  the  normal  distribution  fitting, 
where  in  Fig.  5D  the  assignment  of  channels  G  and  B  with  bands  1 
and  2  resulted  in  the  best  readability  of  vegetation  marks,  while  in 
Fig.  7D  the  use  of  a  single  channel  (band  2)  seemed  to  be  most 
appropriate. 

To  enable  other  colleagues  to  use  the  approaches  discussed  in 
this  paper,  ARCTIS  will  be  released  under  creative  commons  licence 
free  of  charge  via  website  (http://luftbildarchiv.univie.ac.at).  While 
ARCTIS  provides  a  large  bundle  of  functions  for  archaeological 
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Fig.  7.  Data  from  Carnuntum,  Lower  Austria,  acquired  on  June  18  2012.  (A)  Conventional  orthorectified  aerial  image,  acquired  in  the  visible  spectrum.  GSD  of  0.1  m;  enhanced  using 
contrast  limited  adaptive  histogram  equalization  (CLAHE).  (B)  false  colour  composite  created  by  means  of  the  RE1P  algorithm  (R  =  band  1  (wavelength),  G  =  band  2  (slope),  B  =  band 
3(reflectance  value)),  (C)  rate  parameter  b  of  the  fitted  gamma  PDF;  (D)  the  standard  deviation  o  of  the  fitted  normal  PDF.  GSDs  of  (B),  (C)  and  (D):  0.4  m.  Figures  B,  C,  and  D  were 
subject  to  the  same  histogram  stretch  by  means  of  standard  deviation.  (For  interpretation  of  the  references  to  colour  in  this  figure  legend,  the  reader  is  referred  to  the  web  version  of 
this  article.) 


information  extraction,  it  is  not  meant  to  be  a  competitor  for  pro¬ 
fessional  commercial  software  packages.  More  than  anything  else, 
this  toolbox  and  its  residing  algorithms  must  be  considered  a  proto¬ 
typing  and  testing  environment.  It  is  well  known  that  the  MATLAB® 
programming  language  focuses  more  on  the  ease  of  code  devel¬ 
opment  than  on  the  performance  found  in  languages  such  as  C, 
C++  and  Fortran.  This  characteristic  becomes  very  apparent  when 
executing  the  distribution  fitting  approach.  Since  all  operations  are 
executed  on  a  pixel-by-pixel  basis,  the  whole  operation  becomes 
very  time-consuming  in  a  MATLAB®-driven  environment,  which  in 
turn  restrains  it  from  numerous  test  runs  on  large  datasets. 
Although  ARCT1S  was  initially  developed  to  mine  archaeological  AIS 
data  cubes,  its  current  state  prevents  it  from  processing  big  data¬ 
sets.  In  near  future,  when  several  AIS  datasets  have  been  processed 
and  the  outcome  of  specific  processing  workflows  was  validated  by 
different  practitioners,  those  specific  algorithms  (or  combinations 
of  them)  that  proved  most  useful  will  be  implemented  in  C  or  C++ 
for  faster  execution,  while  the  MATLAB®  GUI  can  still  be  used  to  call 
these  external  programs. 

Finally,  archaeologists  also  need  to  find  better  ways  to  quanti¬ 
fying  the  gain  of  image  processing  algorithms.  To  date,  the  common 
practise  in  comparing  different  versions  of  one  image  (or  different 
datasets)  is  by  subjective,  qualitative  means.  Science  and  the 
communication  of  results  often  ask  for  more  objective,  quantitative 
measures  of  possible  gains.  In  addition  to  the  development  of  new 
AIS  processing  algorithms  and  contrast  enhancement  techniques 
for  all  optical  airborne  imagery,  the  authors  are  currently  working 
on  new  approaches  to  more  robustly  express  the  added  value  of 


these  methods.  The  final  aim  of  these  attempts  is  to  maximise  the 
information  extraction  of  remotely  sensed  archaeological  footage, 
but  also  limit  the  AIS  processing  toolbox  to  a  more  confined  set  of 
approaches  that  generally  deliver  the  biggest  gains  for  extracting 
and  interpreting  archaeological  features. 

7.  Conclusion  and  prospect 

This  study  has  investigated  new  approaches  for  analysing 
archaeological  airborne  imaging  spectroscopy  data:  the  red  edge 
inflection  point  (REIP)  and  distribution  fittings.  REIP  extracts  infor¬ 
mation  about  plant  vigour  from  the  steep  sloping  spectral  reflec¬ 
tance  curve  around  700  nrn  (transition  from  the  visible  far-Red  to 
NIR  spectrum).  Finding  the  location  of  the  highest  gradient  in  the 
Red  Edge  spectral  profile  for  each  pixel,  the  algorithm  generates  a 
three  band  image  (wavelength,  slope,  and  reflectance).  The  distri¬ 
bution  fitting  will  for  each  pixel  fit  a  given  PDF  to  the  frequency 
histogram  of  its  spectral  signature  and  report  the  fitted  parameters. 

Both  algorithms  were  programmed  into  a  specific  free  of  charge 
MATLAB®  toolbox  called  ARCTIS,  which  was  created  to  test 
currently  available  AIS  processing  practices  as  well  as  validate  the 
value  of  completely  new  information  extraction  techniques.  Two 
different  high-resolution  datasets,  collected  over  the  same  area 
(the  Roman  town  of  Carnuntum,  Austria)  within  two  years  showing 
different  ground  conditions  were  used  to  test  the  algorithms.  The 
imagery  under  study  has  a  ground-sampling  distance  of  40  cm, 
which  allowed  extraction  of  archaeologically  relevant  features  in 
detail. 
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Fig.  8.  Mapping  of  visualizations  from  Fig.  7:  (A)  orthophotograph,  (B)  RE1P,  (C)  gamma  distribution  fitting,  (D)  normal  distribution  fitting. 


The  examples  given  here  prove  that,  when  compared  to 
simultaneously  acquired  broadband  visible  imagery  at  a  spatial 
resolution  of  10  cm,  distribution  fitting  and  the  REIP  visualisations 
can  increase  the  contrast  between  vigorous  and  stressed  plants. 
That  is,  however,  when  the  data  are  acquired  at  a  high  spatial 
resolution  and  at  the  correct  moment  in  the  plant's  phenological 
cycle.  Once  the  vegetation  marks  are  also  visibly  very  discernable, 
the  specific  advantages  of  both  methods  are  seriously  reduced  and 
the  vertical  photographs  show  most  of  the  archaeological  infor¬ 
mation  in  very  high  detail. 
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